
R version 4.3.1 (2023-06-16) -- "Beagle Scouts"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: aarch64-apple-darwin20 (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(sna)
Loading required package: statnet.common

Attaching package: ‘statnet.common’

The following objects are masked from ‘package:base’:

    attr, order

Loading required package: network

‘network’ 1.18.1 (2023-01-24), part of the Statnet Project
* ‘news(package="network")’ for changes since last version
* ‘citation("network")’ for citation information
* ‘https://statnet.org’ for help, support, and other information

sna: Tools for Social Network Analysis
Version 2.7-1 created on 2023-01-24.
copyright (c) 2005, Carter T. Butts, University of California-Irvine
 For citation information, type citation("sna").
 Type help(package="sna") to get started.

> # Script to run QAP
> library(sna)
> library(doParallel)
Loading required package: foreach
Loading required package: iterators
Loading required package: parallel
> library(doRNG)
Loading required package: rngtools
> 
> # Load the followers adjacency matrix
> y <- readRDS("processed_data/followers_adjacencyMatrix.Rds")
> # Load the predicting matrix 
> load("processed_data/QAP_predicting_matrices.RData")
> 
> diag(y) <- NA
> 
> y_v <- c(y)
> 
> 
> Var.names <-c("State_Similarity",
+               "Party_Similarity",
+               "Chamber_Similarity",
+               "Gender_Similarity",
+               "Race_Similarity",
+               "Difference_in_Legislatures_Profeshionalism",
+               "Dem_Sender_Effect",
+               "Rep_Sender_Effect",
+               "House_Sender_Effect",
+               "Female_Sender_Effect",
+               "Profesh_Sender_Effect",
+               "Black_Sender_Effect",
+               "Latino_Sender_Effect",
+               "Asian_Sender_Effect",
+               "Mena_Sender_Effect",
+               "Multi_Sender_Effect",
+               "Native_Sender_Effect",
+               "Democrat_Receiver_Effect",
+               "Republican_Receiver_Effect",
+               "House_Receiver_Effect",
+               "Female_Receiver_Effect",
+               "Profesh_Receiver_Effect",
+               "Black_Receiver_Effect",
+               "Latino_Receiver_Effect",
+               "Asian_Receiver_Effect",
+               "Mena_Receiver_Effect",
+               "Multi_Receiver_Effect",
+               "Native_Receiver_Effect",
+               "Same_Party_Same_State",
+               "Same_Chamber_Same_State",
+               "Same_Gender_Same_State",
+               "Same_Race_Same_State",
+               "Contiguous_States")
> 
> for(i in 1:length(Var.names)){
+   xi <- predicting_matrices[i,,]
+   diag(xi) <- NA
+   print(sum(is.na(xi)))
+   if(i == 1){
+     xdf <- cbind(c(xi))
+   }
+   if(i > 1){
+     xdf <- cbind(xdf,cbind(c(xi)))
+   }
+ }
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
[1] 4108
> 
> send_id <- matrix(1:nrow(y),nrow(y),nrow(y))
> rec_id <- t(send_id)
> diag(send_id) <- NA
> diag(rec_id) <- NA
> 
> xdf <- cbind(xdf,c(send_id),c(rec_id))
> 
> xdf <- data.frame(xdf)
> 
> names(xdf) <- c(Var.names,c("send_id","rec_id"))
> 
> xydf <- data.frame(y_v,xdf)
> 
> xynd <- na.omit(xydf)
> 
> covs <- paste(Var.names,collapse="+")
> 
> names(xynd)[1] <- "yij"
> 
> library(fastglm)
Loading required package: bigmemory
> system.time(est0 <- fastglm(as.matrix(xynd[,Var.names]),xynd$yij,family=binomial(),method=3))
   user  system elapsed 
 30.694  55.972 158.096 
> 
> estim <- NULL
> twop <- NULL
> for(i in 1:50){
+   filenm =paste("./qap_power_results/qap_power_results/qap_sim_follow",i,".RData",sep="")
+   load(filenm)
+   estim <- cbind(estim,qap_sim$coefficients)
+   twop <- cbind(twop,qap_sim$pgreqabs)
+ }
> 
> psig <- NULL
> for(i in 1:length(coef(est0))){
+   same_sign <- sign(estim[i+1,]) == sign(coef(est0))[i]
+   sig <- mean((twop[i+1,] <= 0.05)*same_sign)
+   psig <- c(psig,sig)
+ }
> 
> # calculate standardized coefficients
> abs_std_coef <- coef(est0)/est0$se
> abs_std_coef <- abs(abs_std_coef)
> 
> names(coef(est0))[abs_std_coef > 2 & psig < 0.8]
[1] "Gender_Similarity"       "Black_Sender_Effect"    
[3] "Female_Receiver_Effect"  "Profesh_Receiver_Effect"
[5] "Native_Receiver_Effect" 
> 
> library(ggplot2)
> # Basic scatter plot
> df <- data.frame(power = psig,zstat = abs_std_coef,
+                  col = ifelse(abs_std_coef > 2 & psig < 0.8,"red","blue"),
+                  shp = ifelse(abs_std_coef > 2 & psig < 0.8,17,19))
> 
> ggplot(df, aes(x=zstat, y=power)) +
+   geom_point(size=2,shape =df$shp,colour=df$col) + 
+   scale_x_continuous(trans='log2') +
+   geom_hline(yintercept=0.8) +
+   geom_vline(xintercept=2) +
+   theme_bw()
> 
> ggsave("power_results_follow.pdf",width=4,height=3.5)
> underpowered <- names(coef(est0))[abs_std_coef > 2 & psig < 0.8]
> write.csv(underpowered,file="underpowered_follow.csv")
> 
> 
> 
> 
> 
> 
> 
> proc.time()
   user  system elapsed 
105.552 129.243 409.886 
